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ABSTRACT 

In order to explain the significant orbital eccentricity of the short-period transiting 
Neptune-mass planet GJ 436b and at the same time satisfy various observational con- 
straints and anomalies, Ribas, Font-Ribera and Beaulieu have proposed the existence 
of an eccentric low- mass companion planet at the position of the outer 2:1 resonance. 
The authors demonstrate the viability of their proposal using point-mass three-body 
integrations, arguing that as long as the system appears to be dynamically stable, the 
short-term secular variations ought to dominate the long-term dissipative evolution. 
Here we demonstrate that if one includes tidal dissipation, both orbits circularize after 
a few times the circularization timescale of the inner planet. We conclude that with 
or without a nearby companion planet, in or out of the 2:1 resonance, the Q- value of 
G J 436b must be near the upper bound estimate for Neptune if the system is as young 
as 1 Gyr, and an order of magnitude higher if the system is as old as 10 Gyr. We 
show detail of passage through resonance and conclude that even out of resonance, a 
companion planet should still be detectable through transit timing variations. 

Key words: planetary systems - celestial mechanics - methods: analytical - plane- 
tary systems: formation 



1 INTRODUCTION 

GJ 436b w as dis covered in a radial velocity survey by 
20041 ') , and has since been observed in transit by 



Butler et a 



Gillon et al 



20071 ) and lAlonso et all (|2008h as the only tran- 
siting hot Neptune to date. Its system parameters are listed 
in Table [T] together with those for the hypothetical compan- 
ion of Ribas, Font-Ribera & Beaulieu (2008) . Here m and R 
are the mass and radius of the body l|Deming et alll2007l ). a 
is the semimajor axis, e is the orbital eccentricity, to, and 
i are the argument of periastron, longitude of the ascending 
node and inclination respectively, the latter two measured 
relative to the line of sight, M is the mean anomaly, Q is 
the Q-value, k is the Love number and R g is the radius of 
gyration. All angles are in degrees; rou nd brackets r efer t o 
the hypothetical orbital data used by I Ribas et all l|2008l) , 
while square brackets refer to initial data used in the simu- 
lation presented here (see Section f3.ip . Note that the mass 
and radius of GJ 436b are respectively 1.35 and 1.094 times 
those of Neptune. As an extremely short-period system with 
a significant eccentricity, it has attracted considerable at- 
tention because it would appear that such a system should 
long ago have circularized. However, there is one factor in 
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particular which works against the intuition that GJ 436b 
ought to be circularized, and that is that the planet's dis- 
tance from the star in units of its radius is quite large at 
160 compared to, for example, HD 209458b at 76. Since 
the circularization timescale depends on the fifth power of 
this ratio, it turns out to be relatively long if one uses es- 
timates for t he Q-value of Neptune. T he latter has been 
estimated bv lBanfield fc Murray! (Il992|) to be in the range 
1.2 x 10 4 < Q N < 3.3 x 10 s , while iTittemore fc Wisdom! 
(1989) estimate the Q-value of Uranus, whose mass is 0.85 
times that of Neptune, to be less than 3.9 x 10 4 . Using the 
expression 



= e b /e b : 



2 (Qb\ (mb\ fOb\ 5 
42tt \k b ) \mj \R b ) b 



(1) 



for the circularization ti mescale of a synchronous system 
l|Goldreich fc Soterl ll966)F1 where the subscript b refers to 
quantities associated with GJ 436b, m* is the stellar mass 
and and P b = 2.64 d is the orbital period, one obtains the 
range 5. 3 x 10 7 yr < Tdrc < 1.5 x 10 9 yr, where we have 
used the iBanfiel'd fc Murray! l|l992t ) estimates for Neptune 



1 Note that iGoldreich fc Soterl dl966t> use a modified Q-value, 
Q' b , which absorbs the Love number such that Q' b = 3Q b /2k b . 
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Table 1. Observed and hypothetical orbital and structural data for a GJ 436 two-planet system 
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as well as Jupiter's quadrupole Love number (0.34). Be- 
ing estimated to lie in the range 1-10 Gyr, the age of the 
system is quite uncertain ([Torres 200?]). Thus in order for 
the non-zero eccentricity to be simply a result of a circu- 
larization time which is longer than the age of the system, 
this simple analysis suggests that Qt must be greater than 
2.3 x 10 5 (r a9e /Gyr), where T age is the age of t he system. 

In this Letter we focus on the proposal of iRibas et alJ 
(2008) that, given the circularization timescale is consid- 
erably less than the age of the system, the eccentricity is 
sustained by the presence of a low-mass companion planet 
positioned at the location of the outer 2:1 resonance, in an 
orbit which is inclined at around 10 to that of GJ 436b. The 
authors claim to have found a strong peak at 5.2 days in a pe- 
riodogram analysis of the RV data, consistent with a body in 
the 2:1 resonance and with a false alarm probability of 20%. 
This is supported by a fit to the residuals of the two-body fit 
to the inner orbit. Moreover, the authors claim that harbour- 
ing a planet in the 2:1 resonance allowed them to produce 
the observed eccentricity of GJ 436b with a planet whose 
mass is low enough not to have previously been detected in 
the RV data0 Their scenario was further strengthened by 
the fact that it provides a natural explanation for the appar- 
ent change in the inclination to the line of sight of the orbit 
of GJ 436b, if,, which wou ld bring it into tran sit between the 
time of the null result of | Butler et al (2004) and the posi- 
tive detection of iGillon et a ]|(|2007fl. The authors estimated 
that a rate of change of it of around 0.1 yr _1 was needed to 
be consistent with these observations (note that this gener- 
ally includes a combinati on of nutation and p recession of the 
orbital plane). However, lAlonso et a] (|2008m have more re- 
cently reported an upper limit for the current rate of change 
of i b of 0.03 ± 0.05° yr" 1 . Since < di b /dt < 0.03° yr^ 1 for 
about 20% of the precession cycle (see Figure HJc)), this in 
itself does not rule out a Ribas et al-type system. 

In the following section we briefly consider the tidal 
evolution of a single-planet system, putting a more accu- 
rate lower bound on Qt for that case. In Sect i on 3 w e study 
the proposed orbital solution of lRibas et al.l (|2008h using a 
direct integration scheme for three bodies which includes 
perturbing accelerations due to tidal dissipation and spin- 
orbit coupling in both the star and the innermost planet, 
and the post-Newt onian relativistic contr ibution to the po- 
tential of the star (|Mardling fc Linll2002l ). We demonstrate 
that long-term evolution tends towards the doubly-circular 
state, independent of the tidal circularization timescale of 
the outer planet, with the orbits inclined to each other by 



2 Note that we find a maximum value of only 0.1 for et over 
a modulation cycle using their initial conditions, independent of 
initial phases and longitudes. 



GJ 436b: Q b =10 5 



0.6 




0.028 1 — 1 — 1 — 1 — 1 — 1 — 1 — 1 — 1 — 1 — 1 — 1 — 1 — 1 — 1 — 1 — 1 
-0.4 -0.2 0.2 0.4 

time (Gyr) 

Figure 1. Past and future evolution of GJ 436b without a com- 
panion. The system cannot be older than 0.38(Q(,/10 5 ) Gyr. 

a fixed angle. Section 4 discusses stability while Section 5 
presents a conclusion. 



2 GJ 436B WITHOUT A COMPANION 

The simplest explanation for the significant orbital eccen- 
tricity of GJ 436b is that the circularization timescale is 
longer than the age of the system. Figure CD shows the re- 
sults of integrating the secular equations l|Mardling fc Linl 
[2002) f° r a single synchronized planet backwards and for- 
wards in time, with system parameters taken from Table [1] 
and Q b = 10 5 (so lid curves) . The equations used are based 
on the analyses of lHutl (|l98ll ) and lEggleton. Kiseleva fc Hutl 
( 1998) in which a constant time lag is assumed for the 
tidal bulge, and are correct for any eccentricity (in as far 
as the constant time lag assumption is correct). If one as- 
sumes that no forces other than those included in the simu- 
lation were acting in the past, and also that the Q- values of 
the planet and star as well as their radii have not changed, 
these results suggest that the system cannot be older than 
0.38(Q p /10 5 )Gyr. This is in contrast to the results obtained 
if one integrates equations which are correct to second-order 
only in the eccentricity (dashed curves). The latter approx- 
imation is clearly adequate for estimating the current cir- 
cularization timescale, but is entirely inadequate for under- 
standing the past evolution. In particular, the second-order 
analysis does not allow an upper bound estimate of the age 
of the system. While alternative tidal models such as those 
which assume constant lag angles for all tidal components 
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jGoldreich fc Soterl Il966'), or tho se which includ e the dy- 
namical tide for high eccentricity (|Mardling||l995h may pro- 
duce slightly different functional dependence on the eccen- 
tricity were the analysis to be carried out, they are not likely 
to change the estimated age of the system by much. 

Other forces are likely to have been operating early in 
the life of the system, including those due to the rapid rota- 
tion of the star, larger stellar and planetary radii, the pres- 
ence of a protoplanetary d isk, coupling of the stella r and 
planetary magnetic fields l|Laine. Lin fc Donel l2008l l. and 
different Q-values. Such factors will have determined the 
"initial" orbital state of the system, that is, the state of the 
system when such influences became negligible, and it seems 
likely that the eccentricity would not have been much larger 
than its present value if the single-planet scenario is correct. 

The main conclusion one can draw from this simple 
analysis is that if GJ 436b does not have companions capa- 
ble of sustaining the observed eccentricity, the upper bound 
on its age in terms of Q b puts a lower bound on the Q- 
value of GJ 436b of Q b > 3 x 10 (r a9e /Gyr), a value slightly 
higher than that gleaned from equation [T] In particular, it 
is approximately equal to the upper bound for the Q-value 
of Neptune if the system is as young as 1 Gyr, and an order 
of magnitude larger if it is as old as 10 Gyr. 



3 LONG-TERM TIDAL EVOLUTION OF 
TWO-PLANET SYSTEMS 

The second simplest explanation for a non-zero eccentric- 
ity in a system like GJ 436 is the presence of a companion 
planet or star. In general the induced eccentricity of a stable 
system will vary quasi-periodically, with periods of variation 
dominated by the rate of change of the angle between the ap- 
sidal lines, the rate of change of the resonance angle(s) if the 
system is in a resonance, and the rate of change of the argu- 
ment of periastron if the sy stem is significantly non-coplanar 
|Murrav fe Dermotd 120001 ). The amplitude of variation de- 
pends on the ratios of semimajor axes and planet masses 
as well as the initial eccentricities and the angle between 
the apsidal lines, and when relativistic effects are impor- 
tant, on the ratio of the inner planet mass to the stellar 
mass. Formulae for the amplitu de and period o f variation 
of the eccentricities are given in IMardlind (|2007T ) for copla- 
nar non-resonant systems with moderate inner eccentricity 
(e& <, 0.2). The expression for the amplitude is also accu- 
rate for moderately non-coplanar and/or resonant systems, 
while that for the modulation period severely overestimates 
the true value for a system in or near the 2:1 resonance 
because for such a close system, higher order and resonant 
terms should be included. 

A mistake often made is that secular variations of the 
orbital elements persist for the lifetime of the system, even 
when dissipative forces are significant (eg. iManess et al.l 
I2007D . In fact such variations are damped out after a few 
circularization timescales, with the system evolving to a 
pseudo-equilibrium configuration which itself evolves on a 
generally longer timescale. A secular theory for two-planet 
coplan ar systems with dissipation has recently been devel- 
oped l |Mardlinell2007t) . in which it is shown that the ec- 
centricity of the innermost planet together with the angle 
between the apsidal lines, r\, evolve towards a fixed point in 
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Figure 2. Possible values for m c and e c corresponding to e^ e? ' = 
0.15 with <2 C = 0.045 AU. No solutions exist for e c < 0.18, while 
essentially all of those with e c > 0.18 are unstable. 



{e b ,rj) space on a timescale of three times the circulariza- 
tion timescale of the inner planet, after which the system 
evolves to the d oubly circular st ate on timescale given by 
equation (60) in IMardlind (|2007l ). Note that this timescale 
is independent of the Q-value of the outmost planet. When 
resonant terms are important and/or the system is moder- 
ately inclined, the behaviour is only slightly modified (as 
long as the system is stable) . The fixed point or equilibrium 
eccentricity, a quantity which is independent of the initial 
values of eh, phases and longitudes, is given by 



>9) 



(5/4) (a b /a c 



yja b /a c {m b /m c )£ c 1 + -ye'i 



(2) 



where e c 



vT 



and 7 = 4(n i) a„/c) 2 (m*/m c )(a c /a„) 3 , 



with n b the mean motion of the inner planet and c the speed 
of light0 Note, however, that ((2)1 was derived assuming that 
the average value of the outer eccentricity doesn't vary much 
on the circularization timescale. When it does (as happens 
for the hypothetical GJ 436 system because the outer mass 
is so low), the estimate ([2| tends to overestimate the equi- 
librium eccentricity. Thus (J5J) can be regarded as an upper 



bound for e 



3.1 GJ 436b with a resonant companion 

If the circularization timescale of GJ 436b is less than a third 
of the age of the system, it will already have reached and 
evolved past the equilibrium eccentricity. Thus we begin by 
calculating the range of companion masses and eccentricities 
capable of producing an equilibrium eccentricity of 0.15, that 
is, the observed eccentricity of GJ 436b. As in the previous 
Section, this will allow us to put a lower bound on its Q- 
value. Using this value in Q as well as a period ratio of 
2, we can rearrange the equation to write m c as a function 
of e c . The result is plotted in Figure [2] which shows that 
no solutions exist for e c < 0.18. Solutions corresponding to 
e c > 0.18 are actually unstable (see next section) so that in 
fact no systems exist for which the equilibrium eccentricity 
is equal to the observed eccentricity. (In fact, the equilibrium 



3 The quantity 7 is the ratio of the relativistic to the three- 
body contribution to the rate of apsidal advance of the innermost 
planet. 
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eccentricity corres ponding to the hypothetical companion of 
iRibas et al.l (|2008h (with e c = 0.2) is 0.06.) 

One can conclude from this that if GJ 436b does have 
a single nearby low-mass companion in or near the 2:1 res- 
onance, the system cannot yet have evolved to its pseudo- 
equilibrium state. Its circularization timescale must there- 
fore be longer than one third of the age of the system, 
thereby providing the weaker constraint on the Q-value of 
GJ 436b than in the single-planet case that Qb > 7.7 x 
10 4 (t a9 e/Gyr). In fact, we can do better than this, and ar- 
gue that after only one circularization timescale the system 
will be in apsidal libration w ith zero minimu m e& and a max- 
imum which depends on e c |Mardlindl2007l ). From stability 
considerations e c can't be much more than around 0.3 (given 
it needs to be able to move safely through the resonance as 
it tidally evolves; see Figure Since this corresponds to a 
maximum value of e b in a libration cycle of 0.12, we must 
have that T c i rc > T age , putting the same lower bound on Qt 
as in the single-planet case. 

We finish this section by demonstrating the behaviour 
discussed above. Using the data in Table [T] we performed 
a direct integration using the code described in the Intro- 
duction (|Mardling fc Lin! 2002). The star's Q-value is an es- 
timate, while the Love numbers k (twice the apsidal mo- 
tion constant) and radii of gyration R g for the star and 
planet correspond t o n = 3 and n = 1 polytropes respec- 
tively (|Sternelll94ll ). the latter often taken to approximate 
the structure of Jupiter. The spin period of the star was 
taken to be 20 days while the spin of the planet was taken 
to be synchronous with the orbital motion. Both were taken 
to be aligned with the orbit normal. 

The secular analysis in iMardlin g) (|2007|) demonstrates 
that varying the Q-value of the inner-most planet merely 
changes the timescale on which the system evolves towards 
the final state while not affecting the local secular oscil- 
lation period. Since tidal dissipation in the planet domi- 
nates the tidal evolution while the inner orbit's eccentricity 
is non-zero, we took Qi = 0.1 in order to show the detail 
of various stages of evolution, and confirmed that the evo- 
lutionary timescale scales linearly with Qb for Q b = 1 and 
10. We also confirmed that the system remains stable for 
2 x 10 5 yr with Qt — 10 6 , consistent with the results of 
Section \3. 2 1 Figure [3] shows 4000 years of evolution, equiva- 
lent to 4 (Q(,/10 ) Gyr for general Qt, except that Qj/O.l as 
many secular oscillations will have occurred in that time. It 
passes through its present (hypothetical) configuration after 
1600(Qb/0.1) yr. Panel (a) shows the evolution of e b (bot- 
tom curve) and e c , with the system entering resonance at 
around 200(Qb/0.1) yr, and entering apsidal libration (panel 
(b)) when e b (temporarily) hits zero at 1000(Qb/0.1) yr. The 
latter occurs on a timescale of 1 r c i rc after which the system 
evolves towards the quasi-equilibrium phase on a timescale 
of 2 Tcirc- The subsequent approach of both eccentricities to 
zero is well approximated by the secular theory and occurs 
on a timescale of around 6r C i rc in this case (equation (60) 
iMardlinel [2007 V Panel (c) shows the approach to a fixed 
value of around 10.3 of the relative inclination. 

Figure [4] shows detail of the passage through resonance. 
While the amplitude of variation of the eccentricity is only 
slightly enhanced (so that estimates provided by the secular 
theory are reasonably accurate), the amplitude of variation 
of the orbital period ratio P c /Pt suddenly increases as the 
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Figure 3. Approach to the doubly circular state of a resonant 
two-planet GJ 436 system. See text for discussion. 
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Figure 4. Detail of the resonant phase of Figure [3] Note that 
the period and amplitude of these variations is independent of 
Qb- See text for discussion. 



system crosses the separatrix and enters the 2:1 resonance, 
with a width of around 0.035. The corresponding variation of 
the orbital period of GJ 436b is around 15 mins, a significant 
increase on the variation outside the resonance, whose effect 
on the transit timing would be easily measured (note that 
the libration period for this system is around 200 days). For 
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procedure described in iMardlind (|2008T ). Initial conditions 
corresonding to unstable systems are indicated in the figure 
by red dots; stable systems are left blank. A clearly defined 
boundary is evident, indicating that a system with a suffi- 
ciently low initial value of e c would be free to tidally evolve 
through the 2:1 resonance without the danger of instability. 
Other resonances are also evident, for example, the 3:2 and 
the 5:2. 



1 1.5 2 2.5 3 

Figure 5. A stability map of the region surrounding the 2:1 res- 
onance for a range of values of e c . 

coplanar systems in the 2:1 resonance with mt, rn c <C m„, 
this variation is given by (Mardling, in preparation) 

SP b /P b = 2[l + (m 6 /m c )a 2/3 ] _1 Sa/a 

~ 2~ 2/3 {m c /m b )Sa, (3) 

where a = P a /Pi, Sa — a — 2, that is, the "distance" from 
resonance, and the approximation holds for m c <C m^Q 
However, the GJ 436 system passes through the resonance 
in around 300(Q 6 /0.1) yr = 0.3(Q 6 /10 5 ) Gyr so unless Q b is 
at least equal to the upper estimate of Neptune's Q-value, 
we would not expect to see the system in resonance now 
(given that it was deposited into the r esonance around the 
time of formation) . lAlonso et aj (|2008h find no obvious de- 
parture from linear ephem eris and conc l ude t hat the pro- 
posed resonant solution of iRibas et al.l (|2008T ) is unlikely. 
Once the system leaves the resonance, the amplitude of vari- 
ation of Pb is given by © with Sa replaced by the quantity 
Sa — V Sa 2 — Act 2 , where Act is the width of the resonance 
(Mardling, in preparation). SPb reduces quickly as the sys- 
tem crosses the resonance, however, even when the system 
becomes doubly circularized SPb is still significant at around 
one minute. Thus we conclude th at even out of reso nance, 
the companion planet proposed bv lRibas et all (S2008) would 
be detectab le through transit timing variations given current 
accuracies (lAlonso et aj|2008h . 

Panel (c) shows the variation of the inclination of 
GJ 436b to the line of sight, ib. The resonance has very 
little effect on the period and amplitude of variation, with 
the average rate of change of ib equal to approximately 

O 1 

0.06 yr . This compa res to the observatio nal upper bound 
of 0.03 ± 0.05° yr" 1 by lAlonso etaj (|2008T ). 

3.2 Stability 

Figure [5] shows a stability map for a system composed of 
GJ 436b and a companion planet with a mass of 4.7 M^, 
aligned periastra, a mutual inclination of 10 , and with ini- 
tial companion eccentricity e c and period P c — aPb indi- 
cated by the position in the plot. Three-body integrations 
were performed, with stability being determined using the 

4 See also iHolman h Murray! feOOSft and lAeol et all ll2005ri but 
note that our expression is independent of the eccentricities, and 
is an approximation to an expression valid for any masses. 



4 CONCLUSION 

The main conclusion from this study is that with or with- 
out a single nearby companion planet, in or out of the 2:1 
resonance, the Q-value of GJ 436b must be greater than the 
upper bound estimate for Neptune if the age of the system 
is around 1 Gyr, and up to an order of magnitude greater 
for an age of up to 10 Gyr. Passage through resonance of 
a Ribas et al. two-planet system occurs on a timescale of 
0.3(Qb/10 5 ) Gyr, remaining stable throughout and beyond. 
However, even it were now no longer in resonance, the com- 
panion planet should still be detectable through transit tim- 
ing variations. 
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